#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <math.h>
#include <pthread.h>

#define MAX_STRING 100
#define EXP_TABLE_SIZE 1000
#define MAX_EXP 6
#define MAX_SENTENCE_LENGTH 1000
#define MAX_CODE_LENGTH 40

char train_file[MAX_STRING];
char output_file[MAX_STRING];
char save_vocab_file[MAX_STRING];
char read_vocab_file[MAX_STRING];

int min_count = 5;
int negative = 5;
int num_threads = 1;
int hs = 0 ;
int window = 8;
int cbow = 1, binary = 0, debug_mode = 2;

typedef float real;
const int vocab_hash_size = 30000000;
int *vocab_hash ;

long long vocab_max_size = 1000;
long long vocab_size;
long long train_words = 0;
long long file_size = 0;
long long iter = 15;
long long classes = 0;
long long word_count_actual = 0;

// 词向量的维度
long long layer1_size = 200;

clock_t start;


struct vocab_word{
    long long cn;
    int *point;
    char *word;
    char *code;
    char codelen;
};

struct vocab_word *vocab;
const int table_size = 1e8;
int *table;

real *expTable;
real alpha = 0.025;
real starting_alpha;
real sample = 1e-4;
real *syn0, *syn1neg, *syn1;

void InitUnigramTable(){
    int a, i;
    double train_words_pow = 0;
    double d1, power = 0.75;
    table = (int *) malloc(table_size * sizeof(int));
    for (a = 0; a< vocab_size; a++) train_words_pow = pow(vocab[a].cn, power);
    i = 0;
    d1 = pow(vocab[i].cn, power) / train_words_pow;
    for (a = 0; a< table_size; a++){
        table[a] = i;
        if(a / (double)table_size > d1){
            i++;
            d1 += pow(vocab[i].cn , power)/ train_words_pow;
        }
        if (i >= vocab_size) i = vocab_size -1;
    }
}


unsigned long long GetWordHash(char *word){
    unsigned long long a, hash = 0;
    for (a = 0; a < strlen(word); a++) hash = hash * 257 + word[a];
    hash = hash %  vocab_hash_size;
    return hash;
}


int AddWordToVocab(char *word){
    unsigned int hash, length = strlen(word) + 1;
    if(length > MAX_STRING) length = MAX_STRING;

    vocab[vocab_size].word = (char *)calloc(length, sizeof(char));
    strcpy(vocab[vocab_size].word, word);
    vocab[vocab_size].cn = 0;

    vocab_size++;
    if(vocab_size +2 >= vocab_max_size){
        vocab_max_size += 1000;
        vocab = (struct vocab_word *)realloc(vocab, vocab_max_size* sizeof(struct vocab_word));
    }

    hash = GetWordHash(word);
    while(vocab_hash[hash] !=-1) hash = (hash + 1) % vocab_hash_size;
    vocab_hash[hash] = vocab_size - 1;
    return vocab_size - 1;
    // hash

}

void ReadWord(char *word, FILE *fin){
    int a = 0, ch ;
    while (!feof(fin)){
        ch = fgetc(fin);
        // \r continue
        if(ch == 13) continue;
        //
        if((ch == ' ') || (ch == '\t') || (ch == '\n')){
            // 换行特殊处理，每遇到一个\n就放一个句子结束符
            if(a > 0 ){
                if(ch == '\n') ungetc(ch, fin);
                break;
            }
            if(ch == '\n'){
                strcpy(word, (char *)"</s>");
                return ;
            } else continue;
        }
        word[a] = ch;
        a++;
        // 如果某个单词大于MAX_STRING 后面的部分会背舍弃;但是指针还是继续往前移动
        if(a >= MAX_STRING - 1 ) a--;
    }
    word[a] = 0;
}

void ReduceVocab(){
    unsigned int hash;
    int a, b = 0;
    for (a = 0; a < vocab_size; a++) {
        if (vocab[a].cn > min_count){
            vocab[b].cn = vocab[a].cn ;
            vocab[b].word = vocab[a].word;
            b++;
        } else free(vocab[a].word);
    }
    vocab_size = b;
    for (a = 0; a< vocab_hash_size; a++) vocab_hash[a] = -1;
    for (a = 0; a< vocab_size; a++){
        hash = GetWordHash(vocab[a].word);
        while(vocab_hash[hash] != -1) hash = (hash + 1) % vocab_hash_size;
        vocab_hash[hash] = a;
    }
}

int VocabCompare(const void *a, const void *b){
    return ((struct vocab_word *)b)->cn - ((struct vocab_word *)a)->cn;
}

void SortVocab(){
    int a, size;
    unsigned int hash;
    qsort(&vocab[1], vocab_size-1, sizeof(struct vocab_word), VocabCompare);

    for (a =0; a < vocab_hash_size; a++) vocab_hash[a] = -1;
    size = vocab_size;
    train_words = 0;

    for (a = 0; a< size; a++){
        if ((vocab[a].cn < min_count) && (a != 0) ){
            vocab_size--;
            free(vocab[a].word);
        }else{
            hash = GetWordHash(vocab[a].word);
            while (vocab_hash[hash] !=-1) hash = (hash + 1) % vocab_hash_size;
            vocab_hash[hash] = a;
            train_words += vocab[a].cn;
        }
    }
    vocab = (struct vocab_word *)realloc(vocab, (vocab_size +1) * sizeof(struct vocab_word));
    for (a = 0; a < vocab_size; a++){
        vocab[a].code = (char *) calloc(MAX_CODE_LENGTH, sizeof(char));
        vocab[a].point = (int *) calloc(MAX_CODE_LENGTH, sizeof(int));
    }
    printf("train words : %lld \n", train_words);
}

int SearchVocab(char *word){
    unsigned int hash = GetWordHash(word);
    while(1){
        if (vocab_hash[hash] == -1) return -1;
        if( !strcmp(word, vocab[vocab_hash[hash]].word)) return vocab_hash[hash];
        hash = (hash+1) % vocab_hash_size;
    }
}

void LearnVocabFromTrainFile(){
    FILE *fin;
    char word[MAX_STRING];
    int i ;
    fin = fopen(train_file, "rb");
    long long a;
    if(fin==NULL){
        printf("ERROR: training data file not found!\n");
        exit(1);
    }

    for (a = 0; a < vocab_hash_size; a++ ) vocab_hash[a] = -1;

    vocab_size = 0;
    AddWordToVocab((char *)"</s>");

    while (1) {
        ReadWord(word, fin);
        if (feof(fin)) break;
        i = SearchVocab(word);
        if (i == -1){
            a = AddWordToVocab(word);
            vocab[a].cn =1;
            // printf("%s: %s\n",vocab[a].word, word);
        }else vocab[i].cn++;
        //if(vocab_size > vocab_hash_size * 0.7) ReduceVocab();
    }

    SortVocab();
    file_size = ftell(fin);
    fclose(fin);
}


void SaveVocab(){
    long long i;
    FILE *fo = fopen(save_vocab_file, "wb");
    for (i = 0; i < vocab_size; i++) fprintf(fo, "%s %lld\n", vocab[i].word, vocab[i].cn);
    fclose(fo);
}
void ReadVocab(){
    FILE *fin;
    fin = fopen(read_vocab_file , "rb");
    if (fin == NULL){
        printf("Vocabulary file not found\n");
        exit(1);
    }
    vocab_size = 0;
    long long a, i= 0;
    char c;
    char word[MAX_STRING];
    for(a = 0; a< vocab_hash_size ; a++) vocab_hash[a] =-1;
    while(1){
        ReadWord(word, fin);
        if(feof(fin)) break;
        a = AddWordToVocab(word);
        fscanf(fin, "%lld%c", &vocab[a].cn, &c);
        i++;
    }
    SortVocab();
    int j = 0;

    for(j = 0; j< 10; j++){
        printf("%s: %lld\n", vocab[j].word, vocab[j].cn);
    }

    fin = fopen(train_file, "rb");
    if (fin == NULL){
        printf("ERROR: training data file not found \n");
        exit(1);
    }
    fseek(fin, 0, SEEK_END);
    file_size = ftell(fin);
    fclose(fin);
}

void CreateBinaryTree(){
    long long a, b, i, min1i, min2i, pos1,pos2, point[MAX_CODE_LENGTH];
    char code[MAX_CODE_LENGTH];
    long long *count = (long long *) calloc(vocab_size *2 +1, sizeof(long long));
    long long *binary = (long long*) calloc(vocab_size*2 + 1, sizeof(long long));
    long long *parent_node = (long long *) calloc(vocab_size *2 + 1, sizeof(long long));
    for (a = 0; a< vocab_size; a++) count[a] = vocab[a].cn;
    for(a = vocab_size; a< vocab_size *2 ; a++) count[a]= 1e15;
    pos1 = vocab_size -1;
    pos2 = vocab_size;

    // 数组越界会发生什么
    for (a = 0; a < vocab_size - 1; a++){
        if (pos1 >= 0){
            if(count[pos1]< count[pos2]){
                min1i = pos1;
                pos1--;
            }else{
                min1i = pos2;
                pos2++;
            }
        }else{
            min1i = pos2;
            pos2++;
        }
        if (pos1 >= 0){
            if(count[pos1]< count[pos2]){
                min2i = pos1;
                pos1--;
            }else{
                min2i = pos2;
                pos2++;
            }
        }else{
            min2i = pos2;
            pos2++;
        }
        // 二叉树节点数量的关系
       count[vocab_size + a] = count[min1i] + count[min2i];
       parent_node[min1i] = vocab_size + a;
       parent_node[min2i] = vocab_size + a;
       binary[min2i] = 1;
    }
    for(a = 0; a< vocab_size; a++){
        b = a;
        i = 0;
        while(1){
            code[i] = binary[b];
            point[i] = b;
            i ++ ;
            b = parent_node[b];
            if(b == vocab_size*2-2) break;
        }
        vocab[a].codelen = i;
        vocab[a].point[0] = vocab_size - 2;
        for( b =0; b < i; b++){
            vocab[a].code[i-b-1] = code[b];
            vocab[a].point[i-b] = point[b] - vocab_size;
        }
    }

    //printf("freee\n");
    //for(a = 0; a< vocab_size; a++){
    //    printf("%d\n", *vocab[a].code);
    //}
    free(count);
    free(binary);
    free(parent_node);
}

void InitNet(){
    long long a, b;
    unsigned long long next_random = 1;
    a = posix_memalign((void **)&syn0, 128, (long long)vocab_size * layer1_size * sizeof(real));
    if(syn0 == NULL){
        printf("Memory allocation failed \n");exit(1);
    }

    if(hs) {
        a = posix_memalign((void **)&syn1, 128, (long long )vocab_size * layer1_size * sizeof(real));
        if (syn1 == NULL){
            printf("Memory allocation failed\n");
            exit(1);
        }
        for(a = 0; a< vocab_size; a++) for (b = 0; b< layer1_size; b++){
            syn1[a*layer1_size+b] = 0;
        }
    }

    if(negative > 0){
        a = posix_memalign((void **)&syn1neg, 128, (long long)vocab_size *layer1_size * sizeof(real));
        if (syn1neg == NULL) {
            printf("Memmory allocation failed\n");
            exit(1);
        }
        for(a = 0; a<vocab_size; a++) for(b = 0; b< layer1_size; b++){
            syn1neg[a*layer1_size+b] = 0;
        }
    }

    // 随机初始化词向量
    for(a = 0; a<vocab_size; a++) for(b =0; b<layer1_size; b++){
        // todo 25214903917 这个数字是什么意思
        next_random = next_random *(unsigned long long)25214903917 + 11;
        syn0[a*layer1_size+b] = (((next_random & 0xFFFF) / (real)65536) - 0.5)/ layer1_size;
    }
    // create Binary Tree
    CreateBinaryTree();

}

int ReadWordIndex(FILE *fin){
    char word[MAX_STRING] ;
    ReadWord(word, fin);
    if (feof(fin)) return -1;
    return SearchVocab(word);
}

void *TrainModelThread(void *id){
    //printf("thread: %d\n", *(int *)id);
    long long a, b, d, cw;
    long long sentence_position = 0;
    long long last_word, sentence_length = 0;
    long long word_count =0, last_word_count = 0, sen[MAX_SENTENCE_LENGTH+1];
    long long l1, l2, c, target, label, local_iter=iter;
    long long word ;
    unsigned long long next_random = (long long)id;
    real f, g;
    clock_t now;

    real *neu1 = (real *)calloc(layer1_size, sizeof(real));
    real *neu1e = (real *) calloc(layer1_size , sizeof(real));

    FILE *fi = fopen(train_file, "rb");
    int tmp = fseek(fi, file_size/(long long)num_threads * (long long )id, SEEK_SET);
    printf("file_size: %lld\n", file_size);
    printf("num_thread : %lld\n", (long long)num_threads);
    printf("id: %lld\n", *(long long *)id);
    printf("seek position : %lld\n", file_size / (long long)num_threads * *(long long *)id);
    printf("fseek result : %d\n", tmp);

    while(1){
        if(word_count - last_word_count > 10000){
            word_count_actual += word_count - last_word_count;
            last_word_count = word_count;
            if((debug_mode > 1)){
                now = clock();
                printf("%cAlpha: %f Progress %.2f%% Words/threads/sec: %.2fk ", 13, alpha,
                        word_count_actual/(real)(iter*train_words +1) *100,
                        word_count_actual / ((real)(now - start +1)/ (real)CLOCKS_PER_SEC *1000)
                        );
                fflush(stdout);
            }
            alpha = starting_alpha * (1-word_count_actual / (real)(iter *train_words +1));
            if(alpha < starting_alpha * 0.0001) alpha = starting_alpha * 0.0001;
        }
        // todo 更新alpha
        //
        // 读一个句子
        if (sentence_length == 0){
            while(1){
                word = ReadWordIndex(fi);
                // 文件结束
                if(feof(fi)) break;
                if(word == -1) continue;
                word_count++;
                // 句子结束
                if (word == 0)break;
                // 采样 的原理是什么呢
                if(sample >0){
                    real ran = (sqrt (vocab[word].cn / (sample * train_words)) + 1) * (sample * train_words)/ vocab[word].cn;
                    next_random = next_random *(unsigned long long )25214903917 + 11;
                    if(ran < (next_random & 0xFFFF) / (real)65535) continue;
                }
                sen[sentence_length] = word;
                sentence_length++;
                // 超过最大长度
                if(sentence_length >=MAX_SENTENCE_LENGTH) break;
            }
            sentence_position = 0;
        }

        if(feof(fi) || (word_count > train_words / num_threads)){
            word_count_actual += word_count - last_word_count;
            printf("feof: %d\n",feof(fi));
            local_iter--;
            if(local_iter == 0) break;
            word_count  = 0;
            last_word_count = 0;
            sentence_length = 0;
            fseek(fi, file_size / (long long)num_threads * *(long long *)id , SEEK_SET);
            //printf("fseek result :%d\n", tmp);
            continue;
        }

        word = sen[sentence_position];
        if (word == -1) continue;
        // 初始化neu1, neu1e
        for(c = 0; c< layer1_size; c++) neu1[c] =0;
        for(c = 0; c<layer1_size; c++) neu1e[c] = 0;
        next_random = next_random * (unsigned long long )25214903917 + 11;
        b = next_random % window;
        // cbow
        if(cbow){
            cw = 0;
            for (a = b; a < window *2 + 1 -b ; a++) if(a != window){
                c = sentence_position - window + a;
                if (c < 0) continue;
                if(c > sentence_length) continue;
                last_word = sen[c];
                if(last_word == -1) continue;
                for(c = 0; c<layer1_size; c++) neu1[c] += syn0[c + last_word * layer1_size];
                cw ++;
            }
            if(cw){
                // 求均值
                for (c = 0; c< layer1_size; c++) neu1[c] /= cw;
                printf("%d\n",vocab[word].codelen);
                for(d=0; d < vocab[word].codelen; d++){
                    printf("%d", vocab[word].code[d]);
                }
                printf("\n");
                printf("in hs\n");
                if(hs) for (d = 0; d< vocab[word].codelen; d++){
                    // printf("%d", vocab[word].code[d]);
                    f = 0;
                    l2 = vocab[word].point[d] * layer1_size;
                    for (c = 0; c < layer1_size; c++) f += neu1[c] * syn1[c+l2];
                    printf("%f\n", f);
                    if (c <= -MAX_EXP) continue;
                    else if (f >= MAX_EXP ) continue;
                    else f = expTable[(int)((f + MAX_EXP) * (EXP_TABLE_SIZE / MAX_EXP /2))] ;
                    //
                    g = (1 - vocab[word].code[d] - f ) * alpha;
                    // propagate output -> hidden
                    for (c = 0; c < layer1_size; c++) neu1e[c] += g * syn1[c+l2];
                    // hidden -> output
                    for(c = 0; c < layer1_size; c++) syn1[c+l2] += g * neu1[c];
                }
                if(negative > 0) for (d = 0; d < negative + 1; d++){
                    if (d ==0){
                        target = word;
                        label = 1;
                    }else{
                        next_random = next_random * (unsigned long long )25214903917 + 11;
                        // 随机选择的概率分布
                        target = table[(next_random >> 16) % table_size];
                        if (target == 0) target = next_random %(vocab_size -1) + 1;
                        if (target == word) continue;
                        label = 0;
                    }
                    l2 = target * layer1_size;
                    f = 0;
                    for (c = 0; c< layer1_size; c++) f += neu1[c] * syn1neg[c+l2];
                    if (f > MAX_EXP) g = (label -1) * alpha;
                    else if (f < -MAX_EXP) g = (label - 0) *alpha;
                    else g = (label - expTable[(int)((f + MAX_EXP) * (EXP_TABLE_SIZE/MAX_EXP / 2))]) * alpha;
                    for (c = 0; c<layer1_size; c++) neu1e[c] += g*syn1neg[c+l2];
                    // 更新隐层参数
                    for (c = 0; c< layer1_size; c++) syn1neg[c+l2] += g*neu1[c];
                }
                // hidden->in
                for (a= b; a< window*2+1 -b; a++) if (a != window){
                    c = sentence_position-window +a;
                    if(c < 0) continue;
                    if(c >= sentence_length) continue;
                    last_word = sen[c];
                    if(last_word == -1) continue;
                    for(c = 0; c<layer1_size; c++) syn0[c+last_word*layer1_size]+= neu1e[c];
                }
            }
        }else{
            // skip gram
            for(a = b; a < window * 2 + 1 - b; a++ )if(a!=window){
                c = sentence_position - window + a;
                if(c<0) continue;
                if(c>sentence_length) continue;
                last_word = sen[c];
                if(last_word == -1) continue;
                l1 = last_word * layer1_size;
                for(c =0; c<layer1_size; c++) neu1e[c] = 0;

                if(hs){
                    for(d=0; d< vocab[word].codelen; d++){
                        f = 0;
                        l2= vocab[word].point[d]*layer1_size;
                        for(c=0; c<layer1_size;c++) f+= syn0[c+l1] * syn1[c+l2];
                        if(f<=-MAX_EXP) continue;
                        else if(f>=MAX_EXP) continue;
                        else f= expTable[(int)((f+MAX_EXP)*(EXP_TABLE_SIZE/MAX_EXP/2))];
                        g = (1-vocab[word].code[d]-f)*alpha;
                        for(c = 0; c< layer1_size; c++) neu1e[c] += g* syn1[c+l2];
                        for(c=0; c<layer1_size;c++) syn1[c+l2] += g* syn0[c+l1];
                    }
                }
                // 负采样
                if(negative>0)for (d =0; d< negative+1; d++){
                    if (d ==0){
                        target = word;
                        label = 1;
                    }else{
                        next_random = next_random *(unsigned long long )25214903917 + 11;
                        target = table[(next_random >> 16 %table_size)];
                        if(target ==0) target = next_random %(vocab_size-1) +1;
                        if(target == word) continue;
                        label = 0;
                    }
                    l2 = target *layer1_size;
                    f = 0;
                    for(c= 0; c < layer1_size; c++) f+= syn0[c+l1] * syn1neg[c+l2];
                    if(f>MAX_EXP) g= (label-1) *alpha;
                    else if(f<-MAX_EXP) g = (label -0) * alpha;
                    else g = (label -expTable[(int)((f+MAX_EXP) * (EXP_TABLE_SIZE/MAX_EXP/2))])*alpha;
                    for(c =0; c< layer1_size; c++) neu1e[c] += g*syn1neg[c+l2];
                    for(c = 0; c<layer1_size; c++) syn1neg[c+l2] += g*syn0[c+l1];
                }
                for(c = 0; c<layer1_size; c++) syn0[c+l1] += neu1e[c];

            }
        }
        sentence_position++;
        if(sentence_position >= sentence_length){
            sentence_length = 0;
            continue;
        }
    }
    fclose(fi);
    free(neu1);
    free(neu1e);
    pthread_exit(NULL);
}

void TrainModel(){
    long a, b ;
    printf("Starting training using file %s\n", train_file);
    //LearnVocabFromTrainFile();
    //SaveVocab();
    starting_alpha = alpha;
    ReadVocab();
    InitNet();
    InitUnigramTable();
    start = clock();
    pthread_t *pt = (pthread_t *) malloc(num_threads * sizeof(pthread_t));

    for (a = 0; a < num_threads; a++){
        pthread_create(&pt[a], NULL, TrainModelThread, (void *)&a);
    }
    for(a = 0; a< num_threads; a++){
        pthread_join(pt[a], NULL);
    }
    FILE *fo;
    fo = fopen(output_file, "wb");
    if (classes ==0){
        fprintf(fo, "%lld %lld\n", vocab_size, layer1_size);
        for(a = 0; a< vocab_size; a++){
            fprintf(fo, "%s ", vocab[a].word);
            if(binary) for(b =0; b<layer1_size; b++) fwrite(&syn0[a * layer1_size + b], sizeof(real), 1, fo);
            else for(b = 0; b< layer1_size; b++) fprintf(fo, "%lf ", syn0[a * layer1_size + b]);
            fprintf(fo, "\n");
        }
    }
}


int main(int argc, char **argv){
    int i;
    strcpy(train_file,"text8");
    strcpy(output_file, "model");
    strcpy(save_vocab_file, "vocab");
    strcpy(read_vocab_file, "vocab");

    vocab = (struct vocab_word *)calloc(vocab_max_size, sizeof(struct vocab_word));
    vocab_hash = (int *) calloc(vocab_hash_size, sizeof(int));
    expTable  = (real *)malloc((EXP_TABLE_SIZE + 1) * sizeof(real));

    // todo 怎么理解exp table
    for (i=0; i < EXP_TABLE_SIZE; i++){
        expTable[i] = exp((i / (real)EXP_TABLE_SIZE*2-1) *MAX_EXP);
        expTable[i] = expTable[i] /(expTable[i] + 1);
    }

    hs = 1;
    negative = 0;
    printf("%d\n", argc);
    printf("%s\n", train_file);
    TrainModel();
    return 0;
}

